Surface-plasmon polaritons in a lattice of metal cylinders 
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Plasmon modes of a two-dimensional lattice of long conducting circular wires are investigated by 
using an embedding technique to solve Maxwell's equations rigorously. The frequency-dependent 
density of states is calculated for various values of the wave vector and the filling fraction. At 
low filling fractions, collective modes are all found to accumulate at the surface-plasmon frequency 



being the bulk plasmon frequency. As the filling fraction increases, the interference 



between the electromagnetic fields generated by localized surface-plasmon polaritons leads to the 
presence of new resonances, whose frequency strongly depends on the interparticle separation. For 
touching wires, a number of multipole resonances fill the spectral range between dipole resonances, 
as occurs in the case of a three-dimensional packing of metal spheres. 

PACS numbers: 42.25.Bs,42.70.Qs,73.22.Lp 
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I. INTRODUCTION 

Over the last few years, a great deal of attention has 
been devoted to the investigation of metallo-dielectric 
photonic crystals, partially stimulated by the appearance 
of remarkable optical properties associated with the ex- 
citation of surface-plasmon polaritons.-i^^i 4 *^ These are 
self-sustained collective excitations corresponding to fluc- 
tuations of the surface electron density^ 

To understand the properties of surface-plasmon po- 
laritons, and their possible role in the optical properties 
of artificially structured metallo-dielectric metamaterials, 
we need the photonic band structure of lattices of metal- 
lic inclusions with frequency-dependent dielectric func- 
tions. Lattices of metallic rods,— cylinders , 10 d 1 i 12 i 13 and 
sphere o 14 ' 15 ! 16 ' 17 have all been studied up to now. Three- 
dimensional (3D) face-centred cubic lattices composed of 
metal spheres, their plasmon modes, and interaction with 
light were studied by Yannopapas et al.^ in the low 
filling fraction regime, and wide bands of multipole res- 
onances were later found to occur as the filling fraction 
increases Flat bands of localized surface plasmons were 
also found to occur in a two-dimensional (2D) periodic 
arrangement of long metallic cylinders^ The impact of 
localized multipole plasmons on the optical absorption 
and the energy loss of 2D and 3D metallo-dielectric crys- 
tals was addressed in Refs. Il8lll9i l20l. by looking at the 
imaginary part of the so-called effective dielectric func- 
tion of the composite. 

In this paper, we consider the evolution of plasmon 
bands in a 2D square lattice made up of long metal- 
lic cylinders embedded in a homogeneous medium, as 
the cylinder size increases from the non-interacting limit 
(small radius compared to the lattice constant) to touch- 
ing cylinders where interactions play a key role. We shall 
consider electromagnetic modes with the wavevector per- 
pendicular to the cylinders, in which case the modes fall 



into two categories: E- (or s-) polarization, in which the 
electric field is parallel to the cylinder axis, and H- (or 
p-) polarization, in which the electric field lies in the 2D 
plane. We shall focus our attention on the plasmonic 
regime, where the radius of the cylinders is much smaller 
than the wavelength of the interacting electromagnetic 
radiation. 

At long wavelengths, the behaviour for .E-polarization 
is the same as that of a homogeneous metal with a re- 
duced plasmon frequency given by v'/o-'p ; 10,18 where / 
is the filling fraction and oj p the bulk plasmon frequency 
of the metal cylinders. Conversely, for i7-polarization 
surface plasmons can be excited, and due to the interfer- 
ence between their electromagnetic fields their frequen- 
cies vary strongly with the cylinder separation, which 
suggests that their frequency can be tuned by changing 
the separation. Here we shall show that for this polariza- 
tion and in the limit of touching spheres multipole reso- 
nances fill the spectral range between dipole resonances, 
as occurs in the case of 3D lattices of metal spheres. 



II. THEORY 

We consider a square array of infinitely long metal 
cylinders embedded in a homogeneous medium of dielec- 
tric constant eo; the metal of the cylinders is described 
by the free-electron Drude dielectric function 21 



lu(lu + i/r) ' 



(1) 



where r represents the inelastic scattering time of bulk 
plasmons. In this section, we describe our theory, start- 
ing off with general considerations and going on to dis- 
cuss the embedding method which we use to calculate 
the dispersion of the electromagnetic modes. 
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A. General considerations 

In the long- wavelength limit, in which the radius of 
the cylinders is much smaller than the wavelength of the 
interacting electromagnetic radiation, we should be able 
to replace our composite material by an effective homo- 
geneous medium of dielectric function e aS (bj), such that 
electromagnetic modes propagate with frequency u> given 
as a function of wavevector k by 



lo = kc/y/e c!f (u>), 



(2) 



where c is the speed of light. Due to the anisotropy of the 
material, the effective dielectric function e eff (w) depends 
on the polarization of the propagating electromagnetic 
waves. 



1. E -polarization 

In the case of electromagnetic waves with the electric 
field polarized along the cylinders (s-polarization) , we 
know from elementary electromagnetism^ that the elec- 
tric field is continuous across the interface. This suggests 
that in the long-wavelength limit the composite will be- 
have as a homogeneous medium with the dielectric func- 
tion of this effective medium given by the weighted aver- 
age of the dielectric functions of the constituents ; 23 ' 24 



e off H = fe(co) + (1 - /)e . 



(3) 



Here, / is the volume fraction of the cylindrical inclu- 
sions. Taking the Drude dielectric function, Eq. ©, and 
cylinders in vacuum (eo = 1), Eq. ([3]) yields 



e eff (w) = 1 - / 



u)(u> + ijr) 



(4) 



This shows that for this polarization the optical response 
of the composite material is expected to be that of free 
electrons in a homogeneous electron gas, but with the 
reduced plasma frequency \fjuj p . 



2. H -polarization 

For electromagnetic waves with the electric field polar- 
ized normal to the cylinders, the electric field is strongly 
modified by the presence of the interfaces. 

In the case of a two-component isotropic system com- 
posed of identical inclusions of dielectric function e(uS) in 
a host medium of dielectric constant eo, the polarization 
P can be easily obtained as follows 



P = /^E to + (i-/)^_lE ,„. 



4tt 



4tt 



(5) 



where E in and E out represent the average electric field 
inside and outside the inclusions, respectively. Assuming 



that our composite material can be replaced by an effec- 
tive homogeneous medium of dielectric function e eff (<x>), 
one can also write 



1 , 



-in 



E, 



(6) 



where E is the macroscopic electric field averaged over 
the composite: 



E = /E,„ + (l-/)E out . 
Eqs. ©-© yield the following relation: 
(e eff - e ) E = f(e - e ) E lB . 



(7) 



(8) 



In the case of a single 2D circular inclusion embed- 
ded in an otherwise homogeneous medium, an elemen- 
tary analysis^ 2 - shows that the electric field E in in the 
interior of the inclusion is 



E n 



u- 1/2 



E. 



(9) 



where 



u=[l-e/e ] 1 . (10) 
Introducing Eqs. © and HO} into Eq. © yields 



e eff (w) = e 



1 



u - 1/2 



(11) 



This equation, which coincides with Eq. (27) of Ref. Il3l . 
describes for a Drude metal the surface-plasmon mode 
u> s = + eo of a single cylinder. 

The interaction among circular inclusions in a host 
medium can be introduced approximately in the frame- 
work of the Maxwell-Garnett (MG) approximation. 25 
The basic assumption of this approach is that the av- 
erage electric field E in within an inclusion located in a 
system of identical inclusions is related to the average 
field E out in the medium outside as in the case of a single 
isolated (noninteracting) inclusion, thereby only dipole 
interactions being taken into account. Hence, in this ap- 
proach the electric field E in is taken to be of the form 
of Eq. © but with the macroscopic electric field E re- 
placed by the electric field E out outside, which together 
with Eqs. © and © yields-^ 



e off (w) = e 



1-/- 



1 



(12) 



with m = (l-/)/2. In the dilute (/ -> 0) limit, m = 1/2 
and the Maxwell-Garnett Eq. $TZ§ yields Eq. (fTTT) . 

As optical absorption by a composite is dictated by 
the poles of e off , an inspection of Eq. (fT2"]l shows that for 
Drude cylinders in vacuum the MG approximation pre- 
dicts optical absorption to occur at \/(l — f)/2u) v , which 
in the dilute (/ — > 0) limit yields the surface-plasmon fre- 
quency uip/y/2. 



3 



Conversely, the energy loss of swift charged particles 
is dictated by the poles of the effective inverse dielectric 
function e~J(uj). From Eq. (fl2|) . one finds: 



1 + /- 



1 



(13) 



with n = (1 + /)/2. Equation (fl~3)) shows that for Drudc 
cylinders in vacuum the MG approximation predicts the 
excitation of surface plasmons at + /)/2w p , which 
in the dilute (/ — > 0) limit yields, as in the case of optical 
absorption, the surface-plasmon frequency u> p / \/2. 



B. Embedding technique 

In the absence of free charges or currents, Maxwell's 
equations reduce to the eigenvalue equation 



V x V x E = e(r)— E, 



(14) 



where e(r) is the spatially varying dielectric function, and 
we assume a magnetic permeability /i = 1 everywhere. 
To solve this equation, we use the embedding technique 
described in Ref. [26| and employed in Ref. [l7| to investi- 
gate plasmon bands in 3D periodic arrangements of metal 
spheres. In this approach, the metallic cylinders are re- 
placed by an embedding potential over their surfaces, 
and the electromagnetic field in between is expanded in a 
given basis set. With the embedding method one has, of 
course, complete freedom in the choice of this basis set; 
we use a plane- wave basis set, which represents a natu- 
ral choice for periodic systems. Because the plane-wave 
basis is used only in the vacuum between the cylinders, 
convergence is very rapid. 

We consider an electromagnetic wave normally inci- 
dent on the structure, i.e., with the wave vector taken 
in the plane perpendicular to the axes of the cylinders, 
which we take to be the z-direction. With this choice, 
the cylinder problem simplifies considerably compared 
with the spheres, because the vector wave-equation can 
be written in scalar form. In the case of .©-polarization, 
the problem reduces to solving the wave-equation for 
E z . The other case is i7-polarization, with the elec- 
tric field perpendicular to the cylinders, in which case 
we have a scalar wave-equation for H z . As usual when 
solving Maxwell's equations, the boundary conditions on 
the fields are that the surface-parallel components of E 
and H are continuous across the surfaces of the cylinders. 
This is taken care of by the embedding potential. 



1. E -polarization 

For ©-polarization, Eq. (Tj"4")) becomes the scalar wave 
equation, 



V 2 E Z = e (r)^E z 



(15) 



In the embedding method, we derive a variational expres- 
sion for the eigenvalue u> 2 /c 2 in terms of a trial electric 
field £ z defined only in the region between the cylinders, 
region /, 



f T drV£* z ■ V£ z - f s dr s f s dv' s £* z (s 



, ,2 as 

e fj dr£* -£ z +c^ J s dr s f s dv' s £* |§£ z 



(16) 

The contribution from each cylinder is replaced by the 
double integral over its surface, containing the frequency- 
dependent embedding potential E(rg, r^; cj 2 /c 2 ), evalu- 
ated at a trial frequency loq. The terms involving the 
derivatives of the embedding potential originally arise 
from volume integrals through the cylinders, but pro- 
vide a first-order correction so that £ is evaluated at the 
estimated frequency u> rather than the trial frequency u>q. 

The embedding potential is defined in terms of the ex- 
act solution E z of the wave equation inside the cylinder, 
at frequency u>o, which matches on to the trial solution 
£ z over the surface of the cylinder, 



dE z (r s ) 
dn s 



= f dr' s Z(r s y s )£ z (r' s ) 
Js 



(17) 



- it gives the exact normal derivative corresponding to 
the trial surface amplitude. 

To solve Eq. (TTB]), the trial function is expanded in 
terms of basis functions F; , 



£ z {r) = y^ejJj(r), i^-(r) = exp(ik; • r). 



(18) 



For a 2D periodic lattice of cylinders, suitable basis func- 
tions are plane waves, where the wave vector k^ is given 
by k, = k + gi, with k the Bloch wave vector (two- 
dimensional, in the plane perpendicular to the cylinders) 
and gi is a 2D reciprocal lattice vector. Substituting into 
Eq. (flU)) and finding the stationary values with respect 
to the coefficients ej gives the matrix eigenvalue equation 



Ae = -?Be, 
c z 



(19) 



where the A and B matrices are given by 



JUj = jf drVf?-VF 3 - jfdrs Jdv^F* U - cu 2 ^j Fj 

°(20) 



and 

Baa 



< o I drF* ■ Fj + c 



: W S ** F '§ F - (21 » 



The integrals in region /, over the 2D unit cell excluding 
the cylinder, are easy to evaluate, 



f k 2 (A - nr 2 ), i = j 



(22) 
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and 



dvF* 



(A-irr 2 ), i=j 
Z7r ' Iki-kjl > 1 / J 



(23) 



- A is the area of the unit cell, and r is the cylinder 
radius. 

From a multipole expansion of the solution of the scalar 
wave equation in cylindrical harmonics and matching 
onto the trial wave function on the surface of the cylinder, 
Eq. (JT7J) yields the following expression for the matrix el- 
ement of the embedding potential: 



dr s / dr' s F*HFj = —2-ktk, exp im(4n — (f>j) 
s Js 



Jm(kir)J m (kjr j 
J m (nr) 



J' m {nr). (24) 



Here, J m (x) represent Bessel Functions of the first 
kindpi <j>i is the angle of basis wave vector k^, and k 
represents the magnitude of the (complex) wave vector 
of the solution of Maxwell's equations inside the cylinder 
at trial frequency ojq, k — y'iwo/c. For this polariza- 
tion, convergence is rapid with respect to the maximum 
harmonic, 77i max , in the expansion. 

The complex dielectric function of the lossy metal 
cylinders, Eq. ([T]), means that the frequencies of solu- 
tions of Maxwell's equations are broadened, and rather 
than solving the eigenvalue equation, Eq. (p~9|) . we find 
the corresponding Green function T(r,r';A) given in re- 
gion / by 



r(r,r';A) = X;r«(A)F i *(r)i!>(r') ) 



with 



k 



(A ik -\B ik )r kj (\) = 5, 



(25) 



(26) 



As we know the frequency at which the Green function 
is evaluated, the embedding potential is evaluated at this 
frequency, and the frequency-derivative terms in A and 
B cancel out in Eq. (f2"6"j) . 

The Green function is related to the spectral density 
of the electric field, 



(27) 



which multiplied by 5u> gives the electric-field intensity 
at point r in this frequency range. As the Green func- 
tion can be expressed in terms of the eigenmodes and 
eigenfrequencies satisfying Eq. (fTS"]) , 



r(r,r';A) 



wf/e>-X 



we see that the spectral density is given by 
2oj 



n(r,u>) 



7TCT 



r e(r)Imr(r, r; lo a /c 2 + iS). 



(28) 



(29) 



Integrating n through region / gives 



7TC 



Y^boTij j drF? • Fj (30) 



nj is the quantity we calculate and plot. 



2. H -polarization 

For iJ-polarization, we solve the scalar equation for 
the 2-component of the magnetic field, 



V 2 H z =e(r)—H z , 



(31) 



and then the requirement that the surface-parallel com- 
ponents of H and E are continuous reduces to the con- 
tinuity of H z and ^dH z /dn s across the surface of the 
cylinders. The magnetic variational principle, analogous 
to Eq. P^|) . is given by 



i /, drVH* z ■ VH Z ~S S dv s j s dv' s n* z (s- 



,2 0E r ' 



ft dvHt -H z + c 2 J s dr s J s dr' s H* z ^ 

where the magnetic embedding potential now satisfies 
1 dH z (v s ) 



(32) 



e dn s 



= f dv' s ^(r Sl r' s )n z (r' s ). 
Js 



(33) 



The matrix elements for determining the eigenmodes 
and Green function are the same as in the case of E- 
polarization, apart from the matrix element of the em- 
bedding potential, which with the definition appropriate 
to the magnetic case, Eq. (1331) . has an extra factor of 1/e, 



dv s / dv'sF^Fj = 
s Js 



2nrn 



exp im{ 



e(w ) ^ 

v J m 



J m (nr) 



J' m (Kr). (34) 



Surface-plasmon polaritons occur with i?-polarization, 
and the maximum harmonic m max imposes a limit to the 
plasmon modes over each cylinder, without which the 
density of states increases without limit at the planar 
surface-plasmon frequency uj p /\/2. Controlling m max en- 
ables us to study the evolution of the plasmon bands in 
detail. Of course there is a physical cut-off for the plas- 
mons due to Landau damping, though this is much larger 
than the m max which we use. Unlike in scattering theory, 
where a maximum multipole moment is also imposed, 
all higher values of m are included in our embedding, 
but with a different boundary condition. In the case of 
-H-polarization, the boundary condition for m > m max 
is that the surface-normal derivative of H z vanishes on 
the surface of the cylinders, corresponding to the surface- 
parallel component of E vanishing. 



Once again we evaluate the electric spectral density, as 
the plasmons are essentially electrostatic in nature. Us- 
ing Maxwell's equations, the integrated spectral density 
is given by 

ni(u) = — - — ~ y ImTij I dvVF* ■ VFj, (35) 

where the Fj's are the basis functions used to expand the 
magnetic field between the cylinders. 

III. RESULTS AND DISCUSSION 

We have applied the embedding formalism described 
above to study a square array of infinitely long Drude 
cylinders in vacuum. We use the dimensionless re- 
duced frequency Co = uja/2irc and reduced wave vector 
k = ka/27r, where a is the lattice constant. The reduced 
cylinder radius is defined as f = 2ixrja. 

First we provide a comparison of embedding with 
finite-difference time-domain (FDTD) calculations by Ito 
and Sakodaj- 2 in which both space and time are dis- 
cretized. They calculated the band structure for H- 
polarization for Drude cylinders in vacuum with reduced 
plasmon frequency Co p — 1, and reduced radius f = f .885. 
Results for this system calculated using embedding, with 
121 and 241 plane wave basis functions and m max = 6 are 
shown in table 1, compared with the FDTD results, taken 
from figure 5 of Ito and SakodaJ^ We see first of all that 
there is excellent agreement, except for band 2, where 
FDTD gives an energy slightly higher than our (con- 
verged) results. We also note that with only 121 plane 
waves the embedding method gives essentially converged 
results, though the individual plasmon bands within the 
band limits given in the tables do change slightly. 



Band 


121 


241 


FDTD 


1 


0.039 


0.039 


0.036 


2 


0.566 


0.566 


0.596 


Plasmons 


0.638 - 0.711 


0.638 - 0.706 


0.642 - 0.718 


3 


1.088 


1.086 


1.082 


4 


1.161 


1.161 


1.164 


5 


1.184 


1.184 


1.182 


1 


0.301 


0.301 


0.300 


2 


0.474 


0.474 


0.515 


Plasmons 


0.632 - 0.720 


0.632 - 0.720 


0.637 - 0.728 


3 


0.921 


0.921 


0.927 



TABLE I: Convergence test for //-polarised electromagnetic 
waves in a square array of Drude cylinders in vacuum with 
Cj p — 1, r — 1.885. Upper table: frequencies of bands at 
k x — 0.05, k y = with 121 and 241 basis functions, compared 
with FDTD results of Ito and Sakodap^. Lower table: the 
same for k x — 0.5, k y — 



For studying the evolution of the plasmon bands, in 
the rest of this paper we choose the reduced plasmon fre- 
quency Cd p = 0.1, which for Al (lo v ~ 15 eV) corresponds 
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FIG. 1: Density of states ni(co) for i?-polarized electromag- 
netic waves in a square array of Drude cylinders in vacuum, 
with f — 0.5 (/ = 0.02). Different curves correspond to in- 
creasing values of k x , from 0.001 to 0.1 (with k y — 0). 



to a lattice constant of a ~ 83 A. This value of Co p gives 
results which are universally applicable in the plasmonic 
regime, where the radius of the cylinders is smaller than 
the wavelength of the electromagnetic radiation. The re- 
duced lifetime r = 2vrcr/a is taken as 1000, broadening 
peaks and enabling us to plot densities of states easily as 
a function of frequency. 



A. B-polarization 

For /^-polarized electromagnetic waves (with the elec- 
tric field parallel to the cylinder axis), our numerical 
results for the integrated density of states rti(io) and 
the band structure co(k) accurately reproduce the re- 
sults expected from Eqs. ((I])-© and, therefore, the re- 
sults reported in Ref. [l3| for this polarization. Figure Q] 
gives ni(u>) for various values of k x at a filling fraction 
/ = 0.02. 2 §. The corresponding electromagnetic band 
structure u>(k), determined from the peak frequencies, 
is shown in Fig. [5] (solid circles) for three values of the 
filling fraction /. This figure shows that in the plas- 
monic regime under consideration the optical response of 
the composite material is indeed that of free electrons 
in a homogeneous electron gas, but with the reduced 
plasma frequency \ffuJ p (solid lines in Fig. [2]) expected 
from Eqs. 

The 2D periodicity of the composite material intro- 
duces frequency gaps in the band structure, but these 
are remarkably small. For cylinders with a reduced ra- 
dius f = 2, corresponding to / = 0.32, the first band 
gap at the centre of the side of the 2D Brillouin zone is 
0.0035 in reduced units, and for touching cylinders with 
/ = 0.79 it is 0.0028. 




FIG. 2: /5-polarization band structure as a function of k x 
(with k y = 0), for cylinders of radius r — 0.5,2.0, and ir 
(touching cylinders) ; the corresponding filling fractions / are: 
0.02, 0.32, and 0.79, respectively. The solid circles represent 
the band structure that we have obtained from the peaks in 
our calculated density of states ni(ui). The solid lines corre- 
spond to a homogeneous electron gas with the reduced plasma 
frequency \ffu p . The dashed line is the light line, ui = ck. 



FIG. 4: Density of states ni(ui) for //-polarized electromag- 
netic waves in a square array of Drude cylinders in vacuum, 
with r = 1.0 (/ = 0.08), wave vector k = (0.001, 0), and vary- 
ing m max : 1 (thin solid line), 2 (stippled line), 3 (dotted line), 
and 4 (thick solid line). 



sionless Mie plasmons have frequencies given by 
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FIG. 3: Density of states nj(w) for //-polarized electromag- 
netic waves in a square array of Drude cylinders in vacuum, 
with f = 1.0 (/ = 0.08), ra max = 4, k y = 0, and varying k x : 
0.001 (thick solid line), 0.05 (thin solid line), and 0.1 (dotted 
line) . 



B. //-polarization 

For //-polarized electromagnetic waves (with the elec- 
tric field normal to the axes of the cylinders) , our results 
are rather similar to those found for a 3D lattice of metal 
spheres,— with almost dispersionless surface plasmons on 
the cylinders and the light line interacting with one of the 
two dipole plasmons. A significant difference is that in 
the case of an isolated Drude metallic sphere, the disper- 



w,= W5r (36) 

where I is the multipole quantum number, whereas for 
an isolated Drude cylinder all the surface plasmons have 
the frequency ui p / \/2, irrespective of quantum number m. 
This holds for lattices with low filling fraction (/—►()), 
and consequently low interaction between the plasmon 
modes, resulting in a narrower range of plasmon bands 
in the case of cylinders. For larger filling fractions, the in- 
teractions between the cylinders dominate, and produce 
a wider range of bands, similar to the spherical case. 



/. Density of states 

We start off by considering the integrated density of 
states. Figure [3] shows ni(u) at three wave vectors, for 
cylinders with f — 1 corresponding to / = 0.08. The 
cylindrical surface plasmons constitute the dispersionless 
peaks near w p /^/2 = 0.071, with the light line dispers- 
ing through the plasmons and interacting with a dipole 
mode. 

In order to study the plasmons in more detail, we plot 
in Fig. |4] the plasmon peaks for cylinders with f = 1 
(as in Fig. [3]) but now at a small wave vector [k x = 
0.001) and various values of m max : 1, 2, 3, and 4. We see 
that at this rather small filling fraction (/ = 0.08), the 
plasmon structure is very narrow - much narrower than 
in the spherical case»il The effect of increasing m max is 
simply to increase the height of the peak at u) = 0.071, 
or u p /y/2. In this sense, the results do not converge 
as m max increases. But this is physically correct - the 
plasmon peak at u> p /y/2 must increase without limit as 





FIG. 5: Density of states m{u>) for //-polarized electromag- 
netic waves in a square array of Drude cylinders in vacuum, 
with f = 2.62 (/ = 0.55), wave vector k = (0.001,0), and 
varying m max : 1 (thin solid line), 6 (stippled line), 10 (dotted 
line), and 12 (thick solid line). 
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to 



FIG. 6: Density of states ni(u>) for //-polarized electromag- 
netic waves in a square array of touching Drude cylinders in 
vacuum, with r = 7r (/ = 0.79), wave vector k = (0.001, 0), 
and varying m max : 1 (thin solid line), 6 (stippled line), 10 
(dotted line), and 12 (thick solid line). 



the number of allowed modes increases. However, the 
rest of the structure clearly converges, as we see from 
Fig 12 

This behaviour is confirmed for larger cylinder radii, 
when the electrostatic interaction between the cylinders 
broadens the plasmon structure. Figure [5] shows the den- 
sity of states at a wave vector close to zero (k = 0.001, as 
in Fig. 2]), for a cylinder of radius r — 2.62 corresponding 
to a packing fraction of / = 0.55. We see that the wings 
of the plasmon structure still converge with m max , and 
once again increasing m max increases the height of the 
central peak. For touching cylinders (Fig. [5]), the density 
of states broadens to give plasmon structure between a 



FIG. 7: //-polarization band structure as a function of k x 
(with k y = 0, i.e., along the TM direction), for cylinders of 
radius f = 1.0 (/ = 0.08). The solid lines and circles represent 
the band structure (converged in m max ) that we have obtained 
from the peaks in our calculated density of states ni(uj). The 
dotted lines represent the u>(k) curves derived from Eq. @ 
with the MG effective dielectric function e e g(o;) of Eq. (|12[) . 



rather low frequency and the bulk plasmon frequency, 
ujp = 0.1. In this case, we cannot distinguish a central 
peak at the surface plasmon frequency, and increasing 
7n max alters the structure within the entire plasmon range 
- though the plasmon limits converge. This is similar to 
the behaviour we found earlier for touching spheres^ 

The size of basis set we need for the convergence of 
these calculations is rather small, about 350 basis func- 
tions for the highest m max under consideration: m max = 4 
in Figs. [3] and SI and m max = 12 in Figs. [5] and [SI 



2. Photonic band structure 

The peaks in the integrated density of states ni(ui) cor- 
respond to normal modes of the system, and plotting the 
frequencies as a function of wave vector we obtain the 
photonic band structures shown in Figs.EHHl Also plot- 
ted in these figures (dotted lines) are the u>(k) curves de- 
rived from Eq. with the MG effective dielectric func- 
tion e eff (u;) of Eq. (fT2|) . The limit of the lower branch oc- 
curs at — f)/2u) p , where the MG e ef! (ui) has a pole, 
and the starting frequency of the upper branch occurs at 
a/(1 + f)/2u> p , where the MG e~ t }{ui) has a pole^ In the 
dilute (/ — > 0) limit, both the limit of the lower branch 
and the starting frequency of the upper branch occur at 
the planar surface-plasmon frequency co p / \2. 

The band structure in Fig. [7] is for a small cylinder 
radius, f = 1 (corresponding to a filling fraction of 
/ = 0.08), and exhibits an infinite number of disper- 
sionless plasmon modes at the planar surface-plasmon 
frequency lo p /\/2. Light (singly degenerate, as we are 
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0.2 0.3 
k 



FIG. 8: As in Fig. [7] but now for cylinders of radius f = 2.0 
(/ = 0.32). This is converged in m max . 



plasmons. 

At the larger filling fraction, / = 0.32, considered in 
Fig. [51 the dipole modes and the light branches move 
further away from the planar surface-plasmon frequency, 
again following closely the MG curves. The multipolc 
plasmon modes of lowest order are slightly broadened 
around uj p / y2, but again there is remarkably little dis- 
persion. There is a considerable difference between these 
results and those for the lattice of spheres with the same 
radius, which shows a much wider spread of multipole 
frequencies, and much more dispersion. 

With touching cylinders (Fig. [5]), all multipole bands 
spread out and MG breaks down badly. The frequency 
range of the multipole plasmon bands is not quite as wide 
as in the case of touching spheres, and the lower band 
edge seems to converge at u) as 0.026 with m max = 12. 
(In the case of spheres, the lower band edge at w k 0.02 
was still dropping at a maximum multipole value of 12.) 



0.2 0.3 
k 



FIG. 9: As in Fig. [3 but now for touching cylinders of radius 
f = n (f — 0.79). m max = 12, and although the band edges 
have converged, the individual bands have not. 



only dealing with TJ-polarization) interacts with one of 
the dipole plasmons to give the light line (represented by 
solid circles) , which is reproduced almost exactly by MG 
theory^ (dashed lines). The second dipole plasmon gives 
the almost dispersionless band joining the upper branch 
of the light line at zero wave vector. The dipole-active 
light line branches span the multipole plasmon branches 
here, all at w p /v2, quite different from the case of a lat- 
tice of metal spheres, where it is only for larger filling 
fractions that the light branches enclose the multipole 



IV. CONCLUSIONS 

We have shown that as in the case of a lattice of metal- 
lic sphere s 17 ' 26 the embedding method provides a very 
economical method of calculating electromagnetic waves 
in a lattice of metallic cylinders^ which unlike other 
methods allows to describe accurately plasmon modes in 
the whole range of filling fractions from the dilute limit to 
the case of touching wires. For the case of ii-polarization, 
where surface plasmons dominate the normal modes, the 
results are quite similar to those obtained for a lattice 
of metallic spheres, except that the frequency range of 
the multipole modes is narrower. For .E-polarization, the 
system behaves almost exactly like a dilute electron gas, 
with a reduced plasmon frequency and only minute band 
gaps at the Brillouin zone boundaries. 

The method and results presented here are for the case 
of the wave vector lying in the plane perpendicular to 
the cylinders. We are also considering the general case 
of arbitrary wave vector, and this will be the subject of 
a subsequent paper. 
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